#include <iostream>
#include <vector>
#define PI 3.1415
#include <math.h>
using namespace std;


double fac(int j)
{
    if (j==0)
    return 1.0;

    double ans=1.0;
    for(;j>0;j--)
    {
        ans = ans*j;
    }
    return ans;
}

string generate_poly_rat(vector<double> &co, vector<double> xi)
{
    string poly;
    string pi;
    int n = xi.size();
    poly += to_string(co[0]);

    for (int i = 1; i < n; i++)
    {
        if (!(co[i] < 0))
            poly += "+";

        pi += ".*(x-(";
        pi += to_string(xi[i-1]);
        pi += "))";

        poly += to_string(co[i]);
        poly += pi;
    }
    return poly;
}


class Newton
{
private:
    vector<vector<double>> diff_table;
    void generate_table(const vector<double> &x)
    {
        int n = x.size();
        double temp;
        for (int i, j = 1; j < n; j++)
        {
            for (i = n - 1; i >= j; i--)
            {
                temp = (diff_table[i].back() - diff_table[i - 1].back()) / (x[i] - x[i - j]);
                diff_table[i].push_back(temp);
            }
        }
    }

public:
    vector<double> *interpolation(const vector<double> &x, const vector<double> &f)
    {
        diff_table.clear();

        vector<double> *p;
        for (auto i : f)
        {
            p = new vector<double>;
            p->push_back(i);
            diff_table.push_back(*p);
        }

        generate_table(x);

        p = new vector<double>;
        for (auto i : diff_table)
        {
            p->push_back(i.back());
        }

        return p;
    }
};

class chebyshev : private Newton
{
public:
    pair<vector<double> *, vector<double> *> *Cheby_interpolation(double f(double), int n)
    {
        vector<double> *x = new vector<double>;
        vector<double> fi;
        for (int i = 0; i < n; i++)
        {
            x->push_back(cos((2 * i + 1) * PI / (2 * n)));
            fi.push_back(f((x->back())));
        }
        pair<vector<double> *, vector<double> *> *p = new pair<vector<double> *, vector<double> *>;

        p->second = x;
        p->first = interpolation(*x, fi);
        return p;
    }
};

class Hermit
{
private:
    vector<vector<double>> diff_table;
    void generate_table(const vector<double> &x)
    {
        int n=x.size();
        for(int j=1;j<n;j++)
        {
            for(int i=n-1;i>=j;i--)
            {
                if(diff_table[i].size()<=j)
                {
                    diff_table[i].push_back(  (diff_table[i][j-1]-diff_table[i-1][j-1])/(x[i]-x[i-j]) );
                }else{
                    diff_table[i][j]=diff_table[i][j]/fac(j);
                }
            }
        }


    }

public:
    pair<vector<double> *, vector<double> *> *interpolation(const vector<double> &x, const vector<vector<double>> &f_diff)
    {
        diff_table.clear();

        vector<double> *p;
        vector<int> degree;


        for (auto i : f_diff)
        {
            degree.push_back(i.size());


            for (int j = 1; j <= degree.back(); j++)
            {
                i.assign(i.begin(),i.begin()+j);
                diff_table.push_back( i );
            }
        }

        pair<vector<double> *, vector<double> *> *q = new pair<vector<double> *, vector<double> *>;
        q->second= new  vector<double> ;
        for(int i = 0;i<degree.size();i++)
        {
            for(int j=0;j<degree[i];j++)
             q->second->push_back(x[i]);
        }

        generate_table(*(q->second));
        q->first = new vector<double>;
        for(auto i:diff_table)
        {
            q->first->push_back(i.back());
        }
        return q;
    }
};